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Recent results relating to ray dynamics in ocean acoustics are reviewed. Attention is focussed on 
long-range propagation in deep ocean environments. For this class of problems, the ray equations 
may be simplified by making use of a one-way formulation in which the range variable appears 
as the independent (time-like) variable. Topics discussed include integrable and nonintegrable ray 
systems, action-angle variables, nonlinear resonances and the KAM theorem, ray chaos, Lyapunov 
exponents, predictability, nondegeneracy violation, ray intensity statistics, semiclassical breakdown, 
wave chaos, and the connection between ray chaos and mode coupling. The Hamiltonian structure 
of the ray equations plays an important role in all of these topics. 

INTRODUCTION 

The chaotic dynamics of ray trajectories in ocean acoustics have been explored in a number of recent 
publications (p]|-[|17|). The purpose of the present paper is to provide a review of results relating to 
this topic. Our exposition is brief but is intended to be self-contained. We introduce a sequence of 
ray-based simplifications to the mathematical description of underwater sound propagation in order 
to get a more complete and clear understanding of the underlying propagation physics, especially in 
range-dependent environments. We consider these simplifications as the starting point of developing 
a quantitative theory. It is our opinion that even full wave simulations cannot be used effectively 
without some understanding of the material described in this paper. 

To make our discussion more concrete, we focus our attention on long-range propagation in 
deep ocean conditions. Theoretical results are emphasized, but with an eye toward analyzing 
measurements. For this reason considerable attention is paid to results that can be applied in the 
presence of complicated (nonperiodic) range-dependent ocean structure. In a separate paper many 
of the results presented and discussed here will be applied to the analysis of a particular data set. 

In the next section we review important preliminary material. First, we introduce ray-based 
solutions to the Helmholtz equation. We then introduce the one-way form of the ray equations 
and the standard parabolic approximation. Finally, the motion of rays in a range-independent 
environment is discussed. For this class of problems the ray equations are integrable and the 
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ray trajectories are most naturally described using action-angle variables. In anticipation of the 
material that follows, the action-angle formalism is introduced. 

Section II focuses on the behavior of rays in range-dependent environments, i.e., on noninte- 
grable ray systems. The action- angle formalism is used here to introduce nonlinear resonances 
and the KAM theorem. The notion of ray chaos is discussed, as are Lyapunov exponents. Our 
discussion of the (well known) limitations on the predictability of isolated chaotic trajectories is 
complemented by a discussion of the (generally unappreciated) stability of families of chaotic tra- 
jectories. Also in this section we discuss the connection between an important condition, known as 
nondegeneracy, and ray stability. 

In section III, we discuss ray intensity statistics and related topics, including the distribution 
of finite range estimates of Lyapunov exponents. The results presented here were first described in 
the analysis of an idealized underwater acoustic problem, but have since been encountered in other 
applications. Problems associated with the important task of connecting ray intensity statistics to 
finite frequency wavefield intensity statistics are discussed. 

In section IV, we provide a more general, but brief, discussion of 'wave chaos' - the study of 
wave systems that, in the ray limit, exhibit chaotic motion. This topic falls slightly outside the 
bounds of providing a review of ray dynamics, but is too important to omit. Strong results relating 
to this topic are dificult to obtain. Much of our discussion focusses on the question of whether 
semiclassical (ray-based) wavefield representations break down at the so-called Ehrenfest range, 
which scales as ln(/) where / is the appropriately nondimensionalized wave frequency. 

In section V, we describe the connection between ray chaos and mode coupling. This work 
builds on well-known results on ray-mode duality. The results described here provide a promising 
means of attacking the wave chaos problem inasmuch as finite frequency effects are built into the 
modal description of the wavefield. 

In the final section, we briefly discuss two issues. First, we discuss the principal shortcoming of 
our current knowledge - our relatively poor understanding of the wave chaos problem. Second, we 
discuss the manner in which ideas relating to deterministic chaos complement and/or conflict with 
more traditional ideas relating to the study of wave propagation in random media. 



This section provides background material that is necessary to understand the material that is 
presented in the sections that follow. Starting with the Helmholtz equation we introduce the ray 
equations and their one-way form, the standard parabolic approximation, and the action-angle 
description of ray motion in range-independent environments. 

A. Waves and rays 

Fixed-frequency (cw) acoustic wavefields satisfy the Helmholtz equation, 



where a = 2nf is the angular frequency of the wavefield and c(r) is the sound speed. In anticipation 
of our later focus on describing sound fields generated by a point source in environments in which 
focusing in the azimuthal direction is assumed to be negligible, we shall assume propagation in 
two space dimensions with r = (z,r) where z is depth and r is range. The so-called short wave 
approximation can be used when 
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o > Vc, 
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i.e., when the acoustic wavelength 2ir/k, where k = cr/c, is smaller than all length scales that 
characterize variations in c. Under such conditions the solution to ([!]) can be written as a sum of 
terms, each representing a locally plane wave, 

<r; ( r)=^i i (r)e^( r ). (3) 

3 

Substitution of the geometric ansatz (^) into the Helmholtz equation (|l]) gives, after collecting 
terms in descending powers of a, the eikonal equation, 

(VT) 2 = c- 2 , (4) 

and the transport equation, 

V(A 2 VT) = 0. (5) 

For notational simplicity we have dropped the subscript j on T and A in (Q) and (||). The solution 
to d|]) can be reduced to the solution to the ray equations, 

dr &H dp &H 

dr dp ' dr dr ' ' 

and 

dT dr 

— = C = p — -H (7) 

dr dr 

where p = VT is the ray slowness (also referred to as the momentum) vector and 

H(p,r) = i(p 2 - C - 2 (r)) =0. (8) 

The independent (time- like) variable r satisfies dr/dl = c where dl = |dr|, and dTC/dr = 0. The 
sum in (Q) is over all ray paths z[r) that connect the source at (zq, 0) and the receiver at (z, r). The 
ray equations (||-||) are seen to have Hamiltonian form, which allows many well known results to 
be applied to the acoustic problem in the short wave limit. Equations (|||8]) describe the so-called 
optical-mechanical analogy of wave propagation in the short wave limit. 

For guided wave propagation in the direction of increasing r, the variable r can be used as the 
independent (time-like) variable and equations (^| - may be rewritten 

dz _ dH dp _ dH 
dr dp ' dr dz ' 

and 

dT t dz u nn\ 
= L = p—-H (10) 

dr dr 

where p = dT/dz is the z-component of the slowness vector, 



H(p,z,r) = —yc 2 (z,r)—p 2 (11) 

is minus the r-component of the slowness vector, and dH/dr = dH/dr. These are the so-called 
one-way ray equations. All subsequent analysis is based on these equations (or a parabolic approx- 
imation to these equations, as described below), rather than the slightly more general equations (|6| 
- ||). Ray angles are defined by the condition dz/dr = tan ip where <~p is measured relative to the 
horizontal. Using (|9|) and ( |Tl|) this reduces to cp = sin ip. An immediate consequence of equations 
(||), independent of the form of H(p, z, r), is d(dz / dr) / dz + d(dp/ dr) / dp = 0. This is a statement of 
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Liouvilles theorem, expressing the incompressibility of flow in phase space (p, z) . The one-way ray 
equations @ have the same form as the equations describing the motion of a relativistic particle 
with mass c _1 and Hamiltonian H. 

The transport equation (||) can be reduced to a statement of constancy of energy flux in ray 
tubes. In a notation appropriate for use with the one-way ray equations, the solution to the 
transport equation, assuming a point source, for the j'-th eigenray can be written 



Aj(z,r) 



A 0j \q 2 i\- 1/2 



(12) 



The matrix element q^i, defined below, describes the spreading of an infinitesimal ray bundle. At 
any fixed r, one has 




(13) 



where the stability matrix 
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Elements of this matrix evolve according to 



^-Q = KQ 

ar 



where Q at r = is the identity matrix, and 



K 



d 2 H 
dzdv 

d 2 H 

dp 2 dzdp 




(15) 



(16) 



At caustics (721 vanishes and the Maslov index /x advances by one unit. (For waves propagating in 
three space dimensions, advances by two units are possible.) At these points diffractive corrections 
to ( |12"D must be applied. The normalization factor Aqj is chosen in such a way that close to the 
source (^) matches the Green's function for the corresponding wave equation, (|l|) or a parabolic 
approximation to the one-way form of this equation. (Aqj is different for Helmholtz equation and 
parabolic equation rays; this small difference is of no consequence in the discussion that follows.) 



B. The parabolic approximation 

The standard parabolic wave equation is 

i fid 2 



v fr-\2!4d^- coU{z > r) r 



where u(z,r) w exp(ikor)^(z, r), ko = ct/cq, 



(17) 



U(z,r) = — [1 



2cq \ c 2 (z, r) J ' 

and Co is a reference sound speed. The corresponding ray equations are (S) and ([H]) with H(p, 
replaced by 

HpE{p,z,r) = yp 2 + U(z,r). 



(18) 



z,r 



(19) 



4 



Ray angles satisfy cqp = tan (p. The parabolic approximation is valid when ray angles are small 
and deviations of c(z, r) from cq are small. The parabolic wave equation (|j~7|) coincides with the 
Schrodinger equation with r playing the role of time and k^ 1 playing the role of Planck's constant 
h. It is useful to note this analogy because many tools that have been developed to study quantum 
chaos, discussed below, can be applied to the study of solutions of (|l7|) or (||) under conditions in 
which the ray equations, (^) or (g), admit chaotic solutions. In this regard it is noteworthy that 
there is no quantum mechanical counterpart to transient sound fields as these are characterized by 
the simultaneous presence of a continuum of ko values. 



C. Integrable ray systems 

It is well known that when the sound speed is a function of depth only and the ocean boundaries are 
surfaces of constant z, the Helmholtz equation (and the parabolic wave equation) admit separable 
solutions. Under the same conditions the ray equations also admit simple solutions that make use 
of action-angle variables, and are said to be integrable. The action-angle formalism is important 
in the material in much of the remainder of the paper, so the essential results are presented here. 

When the sound speed is a function of depth only, ray trajectories are periodic; the ray equations 
(0) can be transformed, via a canonical transformation, to a new set of ray equations in which the 
Hamiltonian H (I) is a function of the new momentum variable I but is independent of the new 
generalized coordinate 6. In terms of the action-angle variables (/, 6), ray trajectories are described 
by the equations 

dB dH /T dIdH 

T r = W s ^Tr = -m =°- (20) 

The solution to these equations is simply / = constant, 0(r) = tu(I)r + #(0). These equations 
(and their higher dimensional counterparts) describe motion on a torus. The action variable can 
be written (see, e.g., Ref. jL8|) as a function of H (which is constant following each ray trajectory 
in a range- independent environment), 



1 r l r £ ( H ) 

I= — fdzp(H,z) = - dzp(H,z). (21) 

2ir J ir J 2(H) 



The integration is over one cycle of the periodic ray trajectory, and at the turning depths c 1 (z(H)) = 
c (z(H)) = —H. The generating function for the canonical transformation from (p, z) to (1,0) is 

G(z,I) = J" dz'p{H{I),z') (22) 

where the H(I) can be obtained by inverting equation (|2l|), p = dG(z, I)/dz, and 6 = dG(z, I)/dI. 



II. NONINTEGRABLE RAY SYSTEMS 



Realistic ocean structure is range-dependent, and the corresponding ray systems are nonintegrable. 
Thus, from a practical point of view, one is forced to consider nonintegrable ray systems. Also, from 
a more abstract dynamical systems viewpoint, the property of integrability (separable solutions) 
is known to be rather special. It is now customary to consider nonintegrable systems as being 
'generic'. In other words, if an 'arbitrary' set of governing equations of motion were considered, it 
would most likely be nonintegrable. Only under rare circumstances would a system be found to 
be integrable. The other end of the dynamical spectrum is the completely chaotic system which 
is discussed in more detail ahead. Another possibility includes near-integrable dynamics, so-called 
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because perturbation theory on the nearby integrable system would be valid. Finally, there is 
mixed dynamics in which both regular and chaotic motion simultaneously exist in the system that 
are accessed via different initial conditions. Models of ray dynamics in long-range ocean acoustic 
propagation tend to be in the mixed to fully chaotic regimes. 



A. Nonlinear resonances and the KAM theorem 



Consider ray motion in an environment consisting of a range-independent sound speed profile to 
which a small range-dependent perturbation is added. Because of the smallness of the sound speed 
perturbation, the perturbation to 77 may be assumed to be additive, 

H(p,z,r) = H (p,z) +sH 1 (p,z,r). (23) 

For simplicity, first consider the case where H\ is a periodic function of r with wavelength A = 2tt/Q. 
It is well known (see, e.g., Ref. §) that for this class of problems canonical perturbation theory 
fails as nonlinear ray-medium resonances are excited for those rays whose action values Io in the 
unperturbed environment satisfy the condition 

Iuj(Iq) = mVt (24) 

for any pair of integers, I and m. 

A simple analysis (see, for example, Ref. Q) shows that action variables of rays captured into 
the resonance belong to the interval 7 — A7 max < 7 < 7 + A7 max with 



A7 max = 2^/efli/|w'|, (25) 

where uj' = duj(l)/dl at 7 = 7o, and 77i is a measure of the magnitude of 77i. The quantity A7 max 
represents the width of the resonance in terms of the action variable. The width of the resonance 
in terms of spatial frequency can be approximately estimated as 



Auj = \u'\ A7 max /2 = ye77i \uj'\. (26) 

The phenomenon of nonlinear ray medium resonance plays an important role in the emergence 
of ray chaos. If there are at least two nonlinear resonances centered at spatial frequencies uj and 
uj + 5uj, chaotic motion, according to Chirikov's criterion ( |]l9| -[^lll), takes place when the condition 

> 1 (27) 

OUJ 

is satisfied, i.e., when the resonances overlap leading to the stochastic instability of the system. 

One might expect that, even for very small e, all rays are captured into a nearby resonance. This 
turns out not to be the case. According to the KAM theorem (see, e.g., Ref. ||22|), for sufficiently 
small e some of the tori of the unperturbed system are preserved in the perturbed system, albeit 
in a slightly distorted form. A condition for the applicability of the KAM theorem is that the 
nondegeneracy condition uj' ^ must be satisfied. This condition guarantees that resonances are 
isolated provided e is sufficiently small. Nondegeneracy violation will be discussed in more detail 
below. 

Because realistic sound speed structure in the ocean does not have periodic range-dependence 
it is important to consider a larger class of perturbation terms H\(p, z,r). In Ref. [^3[ it is shown 
that the KAM theorem applies to problems for which the perturbation term 77i consists of a 
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superposition of N components, each of which is periodic in range. N is assumed to be finite 
but is otherwise unrestricted. The spatial periods need not be commensurable and there may be 
depth structure associated with each periodic component. Realistic ocean sound speed structure 
- internal- wave- induced perturbations [24], for example - can be described by a model of this 



type. A consequence is that the mixture of chaotic and regular trajectories (discussed below) that 
characterizes ray motion in environments with periodic range-dependence applies to a much larger 
- and more oceanographically realistic - class of problems. 

We note that there is a strong reason to study ray motion in periodic environments even if 
they seem somewhat unrealistic. In such environments phase space structure can be observed by 
plotting pairs of points (p, z) at integer multiples of the wavelength A of the perturbation. Such a 
Poincare map, which is found numerically, is a slice of the ray motion in (p, z, r mod A). In spite of 
the artificial nature of the periodicity, the main features found are expected to be qualitatively just 
like those found with aperiodic environments given that the KAM theorem is equally valid whether 
the environment has the property c(z, r + A) = c(z, r) for some A or not. 

For some special problems ray dynamics and associated phase space structures can be studied 
using an even simpler technique which eliminates the need to numerically trace rays. An example 
is described in Ref. ||. Here, ray motion in a bilinear model (constant sound speed gradient 
above and below the sound channel axis) of the deep ocean sound channel was considered. It was 
shown that when the upper ocean sound speed gradient oscillates periodically in range, successive 
(separated by one ray cycle) iterates of axial ray angle and range satisfy 

4>n+i = 4>n + e[sin p n + sin(p n + (p n + e sin p n )\ , 

Pn+l = p n + (f>n + £ sin p n + 70 n+ i . (28) 

Here e is the dimensionless perturbation strength, 7 is the ratio of the average upper ocean sound 
speed gradient g to the fixed lower ocean sound speed gradient, <p n = (4tt / 'g\)(p n and p n = (27r/A)r n . 
These equations define an area-preserving mapping, d(4> n+ i, p n+ \) / d(<p n , p n ) = 1; this condition is a 
discrete analog of Liouville's theorem. Because of their relative simplicity, area-preserving mappings 
are widely used to study properties of nonintegrable Hamiltonian systems. Some of these properties 
will now be described. 



B. Ray chaos 

Figure 1(a) shows iterates of equation ([28|) for many sets of ray initial conditions in an environment 
with moderately strong (e = 0.15) range-dependence. Phase space is seen to consist of a mixture 
of regular and chaotic regions, often described as regular islands in a chaotic sea. For much 
smaller values of e, chaotic regions occupy only thin isolated bands of phase space. Each such 
thin chaotic band is associated with an isolated resonance. As e is increased, the widths of the 
resonances increase and nearby resonances overlap leading to an intricate mixture of regular and 
chaotic regions, as seen in Fig. 1(a). Behavior of this type is typical of systems that are constrained 
by KAM theory. Figure 1(b) shows a plot of range vs. launch angle after 250 ray cycles in the 
same environment, and using rays emanating from the same fixed point that was used to produce 
Fig. 1(a). The ray initial conditions used in both Figs. 1(a) and 1(b) fall on a horizontal line (at 
r = A/2) through the middle of Fig. 1(a). Note that the islands that intersect this line in Fig. 1(a) 
are readily identifiable in Fig. 1(b). This observation is significant because plots like Fig. 1(b) - 
ray position vs. some continuous ray label - can be constructed in environments with nonperiodic 
range-dependence, providing a simple means of identifying island-like structures. At the boundaries 
of chaotic and regular domains - in Fig. 1(a), for example - are usually invisible cantori. These 
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angle (deg) launch angle (deg) 

Figure 1: Numerical simulations based on the area-preserving mapping (|28| ) with g = 1/(30 km), 
A = 10 km, 7 = 4 and e = 0.15. In both plots the ray initial conditions correspond to those 
of an axial point source at r = A/2. Left panel: 500 iterates of the mapping for 50 rays whose 
launch angles are uniformly distributed between 8 and 9 degrees. Right panel: range after 250 ray 
cycles (each corresponding to one iteration of the mapping) for 10,000 rays whose launch angles 
are uniformly distributed between 8 and 9 degrees. In most regions this sampling interval is too 
large to resolve what should be a smooth function. 



are Cantor-type invariant sets with fractal structure containing an infinity of holes. Cantori act as 
partial barriers that inhibit the diffusion of rays. Boundaries of chaotic regions contain small island 
chains around which the density of points is very high. A common phenomenon is "stickiness" 
of island boundaries; after wandering in an apparently random fashion in phase space, a chaotic 
trajectory may approach a stable island, and stick to its border for some (possibly long) time, 



during which it exhibits almost regular behavior [25]. The presence of regular islands in phase 
space alters the dispersion characteristics of the trajectories that lie in the surrounding chaotic sea 
[p5| . Details depend on the structure of phase space but the phenomenon of anomalous diffusion 
(mean square displacements of trajectories in phase space obeying scaling laws different than that 
of a traditional random walk) seems to be generic. 

An important property of chaotic trajectories is that they exhibit extreme sensitivity, charac- 
terized by a positive Lyapunov exponent, 

lim (\ lim T>(r)\ 
VL = ™ [r P(0)-0 ln W))- (29) 

Here D(r) is a measure of the separation between rays at range r. A closely related quantity is 
the Kolmogorov-Sinai entropy h^s- Loosely speaking, hxs is a measure of information increase 
following a trajectory; a readable discussion of this topic can be found in Ref. J26[ . For bounded 
dynamical systems h^s ~ V L holds; in open diffusive systems, their difference is proportional to a 
diffusion constant [27]. A consequence of the extreme sensitivity of chaotic rays is that the number 



of eigenrays connecting a fixed source and receiver grows exponentially, like exp(/i^5r), on average 
in range |l], [||, O. Also, the magnitudes of the variational quantities qij (see equation |l4|) can 



8 




Figure 2: Left panel: segment of a Lagrangian manifold for a fan of rays with z(0) = —1.100 
km, 8.575° < 8(0) < 8.625° at r = 500 km in the background environment shown in Fig. 3 
with an internal-wave-induced sound speed perturbation superimposed. Right panel: segment of a 
Lagrangian manifold for a fan of rays with z(0) = —1.101 km, 8.575° < 9(0) < 8.625° at r = 500 
km in the same environment. The endpoints of both Lagrangian manifolds are marked with small 
open circles. In both panels, portions of the manifold segments, consisting of long thin tendrils, 
extend beyond the plot boundaries. The dashed curve in the right panel is a portion of a surface 
/ = constant. 



grow exponentially, on average, in range. Reference |17J contains a detailed discussion of this topic. 
A consequence of this exponential growth is that the amplitudes of chaotic rays (proportional to 
\dz(r)/dp(0)\- 1 / 2 for © or \dp n / d^ 1 ' 2 for (|28|), assuming a point source) decay exponentially, 
on average, in range. Note, however, that for moderate to large r or n (measured in units of a 
typical value of v^ 1 ), plots of dz(r) / dp(0) vs. p(0) (see Fig. 5) or dp n /d(j)o vs. 4>o (see Fig. 1) have 
fractal-like structure when both chaotic and regular trajectories are present. 

Another consequence of extreme sensitivity of chaotic rays is that deterministic predictions 
using finite precision numerical arithmetic is limited to ranges less than some threshhold (which 
is proportional to N^/vl where N), is the number of bits used to specify the mantissa of floating 
point numbers). This limitation on one's ability to make deterministic predictions of chaotic ray 
trajectories at long range is tempered by two related factors. First, the shadowing lemma [26| 
guarantees that, for a large class of problems, numerically computed chaotic trajectories at long 
range correspond to trajectories of rays whose initial conditions are close to the specified values but 
are generally unknown. Second, it is easy to verify that statistical properties of discretely sampled 
distributions of chaotic rays (occupying, for example, a small but finite area in phase space at 
r = 0) evolve in a way that does not exhibit extreme sensitivity. 

A somewhat stronger form of stability of distributions of chaotic rays is illustrated in Fig 2. 
In phase space an aperture-limited compact source is represented as a line segment at constant 
depth z = zq bounded by limiting values of po> this is an example of a Lagrangian manifold. Each 
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1.48 1.50 1.52 1.54 

c (km/sec) 



Figure 3: Background sound speed profile used to produce Figs. 2, 4 and 5. 



point on such a manifold evolves in r according to the ray equations @ . As a Lagrangian manifold 
evolves, it gets stretched and folded, but does so without breaking or intersecting itself (owing to 
phase space area preservation). Fig. 2 illustrates one aspect of a phenomenon that can be termed 
manifold stability. In this figure two Lagrangian manifolds with slightly different initial conditions 
are plotted in phase space at a fixed range, r = 500 km. Each manifold has initial conditions 
z = zq (a constant), 8.575° < 8q < 8.625°. For one of the manifolds zq = 1.100 km; for the other 
zq = 1.101 km. Both manifolds evolved in the same environment, which is described below. In this 
environment almost all trajectories evolve chaotically, which leads to exponential growth in range 
of the length of each manifold. The combination of chaotic ray motion and the small difference 
in initial manifold depth might lead one to expect that the two manifolds should evolve very 
differently. Figure 2 shows clearly, however, that they have not. This can be explained by noting 
that phase space area is preserved. Because of this constraint, the exponential stretching of each 
manifold is balanced by an exponential contraction of phase space in the transverse direction. This 
contraction causes the two manifolds to get squeezed closer to one another. Note, however, that, 
in general, points on the two manifolds with the same value of 9q do not lie close to one another. 
Loosely speaking, the two manifolds have slid relative to one another while being stretched, folded 
and squeezed toward one another. Trajectories with nearby initial conditions will also be squeezed 
toward the same locus of points. Thus, although individual trajectories exhibit extreme sensitivity 
under chaotic conditions, the associated exponential contraction of phase space elements imparts 
a surprisingly strong form of stability on continuous distributions of trajectories. The impact of 
manifold stability on wavefield evolution and stability is considered in Ref. [28]. 

The environment used to produce Fig. 2 is also used in subsequent numerical work. It consists of 
a range-independent background profile, shown in Fig. 3, on which an internal-wave-induced sound 
speed perturbation field is superimposed. The background profile is a Munk |3l[] profile modified 
in the upper ocean, c(z) = co(l + e(exp?7 — 77 — 1)) + c u (z), with cq = 1.49 km/s, e = 0.0057, 
rj = (z — z a )/B, B = 1 km, z a = —1.1 km, and c u (z) = 5sin 2 (-7r(z — z a )/z a ) for z > z a with 



10 



0.018 km/s. The internal-wave-induced sound speed perturbation was computed using Eq. 19 



of Ref. [24] with y = t = 0, i.e., a frozen vertical slice of the internal wave field was assumed. An 
exponential buoyancy frequency N{z) = Nq exp(z/B) (note that depths z are negative with z = at 
the sea surface) with Nq = 6 cycles/hour was used. The dimensionless strength parameters E and 
fi were taken to be 6.3 x 10~ 5 and 17.3, respectively. Numerically, a 2 14 point FFT was used with 
Ak x = 27r/1638.4 km, k x . max = 2-k/I km and j m ax = 30. It should be noted that this perturbation 
field is highly structured and fairly realistically describes typical deep ocean environments. Also, 
the assumed background profile is similar to profiles found in much of the North Atlantic Ocean. 

We have seen that the generic mixture of chaotic and nonchaotic ray trajectories, the limited 
predictability of chaotic trajectories, and the strong constraining influence of Liouville's theorem 
lead to very sublte ray dynamics that blur the distinction between traditional deterministic and 
stochastic viewpoints. These and related ideas will be further explored in sections that follow. 

C. Nondegeneracy violation 

Identifying conditions under which the KAM theorem does not apply is important for two reasons. 
First, inapplicability of KAM theory allows rays to behave chaotically even for infinitesimally small 
perturbations to the background sound speed structure. Second, the character of this chaotic 
motion is expected to differ from the chaotic motion of rays constrained by KAM theory. This 
is because in systems constrained by KAM theory the presence of islands - even small islands 
- constrains the motion of nearby chaotic trajectories. When KAM theory does not apply, the 
constraining influence of islands is lost. In this subsection we briefly discuss a cause for KAM 
theory not to apply - violation of the nondegeneracy condition. We expect that nondegeneracy 
violation is important in a wide variety of underwater acoustic environments, but the topic has not 
yet been systematically explored. 



Ref. [21] makes a distinction between intrinsic and accidental degeneracy. In the former case 
dw/dl vanishes for all trajectories. Parabolic equation rays in the environment (cq/c(z)) 2 = 1 — 
{{z — zo)/h) 2 and Helmholtz equation rays in the environment c{z)/cq = cosh((z — zo)/h) are 
intrinsically degenerate. These special problems are not realistic oceanographically. Also, it is 
worth noting that perturbations to intrinsically degenerate systems often break up the degeneracy, 
leading to phase space structure that is essentially the same as that seen in nondegenerate systems 
[ pl| . Accidental degeneracies, on the other hand, occur when diu/dl has isolated zeros - when duj/dl 
is plotted as a function of launch angle for a point source, for example. Accidental degeneracies 
are structurally stable and are fairly common in realistic models of background ocean structure. 
Thus we expect that accidental degeneracies are important in ocean acoustics. An example of an 
accidental degeneracy is briefly discussed below. 

The foregoing discussion can be quantified somewhat. Recall that two conditions for the validity 
of KAM theory are: 1) the dimensionless perturbation strength in equation (23) e is small; and 2) 
the nondegeneracy condition, 

(30, 

is satisfied in the background environment, which is assumed to be range- independent. Zaslavsky 
[29 has argued that these conditions may be related; he argued that KAM theory is valid provided 

e<C|a| 7 (31) 

(for small e and a) for some system-dependent C and 7 > where a is the dimensionless stability 
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Figure 4: Stability parameter a vs. axial ray angle in the environment shown in Fig. 3. 



parameter 

I du> , . 

a = ~—. (32) 
ijj dl 

The term 'weak chaos' has been used |3C, 0] to describe chaos that is caused by an infinitesimal 



perturbation to an integrable system. When a vanishes, equation (|3lJ) will be violated for arbitrarily 
small nonzero e and weak chaos is expected to be present. In Ref. || it was demonstrated that 
rays in a small band surrounding the surface-grazing ray exhibit weak chaos; for the surface grazing 
ray du/dl - and hence also a - is undefined. 

The stability parameter a can be computed from more familiar ray quantities. The angular 
frequency u = 2ir/R where R is the ray cycle distance. For Helmholtz equation rays 

R(p r ) = 2 Pr J i (c _ 2(z) _^ )1/2 (33) 

where p r is the r-component of the ray slowness vector, which is constant following a ray trajectory 
in a range-independent environment. The upper turning depth of a ray z(p r ) satisfies c(z(p r )) = 
p^ 1 , and similarly for the lower turning depth z(p r ). Also, cp r = costp. The action can be written 

I(Pr) = ~ f dz{c- 2 {z)-p 2 r fl 2 . (34) 

7T Jz 

Because 2irdl/dp r = —R(p r ), a = (2ttI j 1 R?)dR/dp T . This expression for a is convenient to evaluate 
numerically. 

Figure 4 shows a as a function of axial ray angle in the environment shown in Fig. 3. It is seen 
that a has two zero crossings for rays with axial angles near 12°. Figure 5 shows a plot of ray depth 
vs. launch angle for an axial source at a range of 1000 km in this environment with an internal- 
wave-induced sound speed perturbation superimposed. (The internal-wave-induced perturbation 
used to produce Fig. 5 is identical to that used to produce Fig. 2 apart from a factor of two 
reduction in the strength parameter E.) In Fig. 5 it is seen that rays in the 10° to 12° band are 
evidently much less stable than those outside this band. Outside of this band there is evidence 
of apparently nonchaotic island-like structures (compare to Fig. 1); inside this band there is no 
indication - at the resolution shown, at least - of the presence of such structures. We attribute this 
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Figure 5: Ray depth vs. launch angle for an axial source at a range of 1000 km in the background 
environment shown in Fig. 3 with an internal-wave-induced sound speed perturbation superim- 
posed. In most regions the sampling interval (A6q = 0.0005°) is too large to resolve what should 
be a smooth function. 



difference to rays in the more stable regions being constrained by KAM theory. (The slight upward 
shift of the unstable band of axial ray angles in Fig. 4 relative to Fig. 5 is not unexpected. In the 
presence of internal waves in the deep ocean, the upper turning depths of rays whose axial angles 
are more than a few degrees tend to move slightly toward the surface - or equivalently, axial ray 
angles increase slightly - with increasing range. This is because, following such a ray, scattering 
is strongest near the ray's upper turning depth; any perturbation to the ray angle at the upper 
turning depth will make the ray steeper.) 

Because accidental degeneracies are not uncommon in realistic models of (background) ocean 
structure, we expect the qualitative behavior exhibited in Figs. 4 and 5 to be fairly common 
in ocean acoustics. The occurrence of accidental degeneracies may be the cause of the finding of 
Simmen et al. [14] that, for identical sound speed perturbation fields, characteristics of ray stability 
vs. launch angle curves depend strongly on the background sound speed profile. The connection 
between nondegeneracy violation and ray stability needs to be further investigated. 



III. RAY INTENSITY STATISTICS 



Ray intensity distributions are discussed in this section. There are two principal reasons for investi- 
gating this topic. First, ray intensity distributions are a useful diagnostic tool to study and quantify 
ray dynamics, especially in environments with nonperiodic range-dependence where Poincare maps 
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cannot be constructed. Second, one expects - on the basis of the local plane wave expansion (|j) 
- that a useful starting point for understanding wavefield intensity statistics is understanding ray 
intensity statistics. We emphasize that gaining an understanding of ray intensity statistics is only 
the first step in this process, as diffraction and interference effects must be accounted for in making 
the transition to wavefield intensity statistics. These complications will be discussed in more detail 
below. 

We shall confine our attention to a discussion of ray intensity statistics in a simple idealized 
ocean model. A more detailed account of the results presented here can be found in [17]. Some 



related results for the same idealized problem are presented in [16]. The model consists of an 
unbounded homogeneous background on which an isotropic, single scale (described by a Gaussian 
spectrum) random perturbation is superimposed. It should not be expected that all of the prop- 
erties of ray intensity distributions in our idealized problem carry over to long-range deep ocean 
propagation conditions, which are characterized by the presence of a background sound channel on 
which strongly inhomogeneous and isotropic internal-wave-induced sound speed fluctuations with 
a power-law spectrum are superimposed. In spite of the idealized nature of the problem considered 
here, there is ample justification for studying this problem. First, this simple problem has surpris- 
ingly rich structure that must be understood before more complicated problems can be successfully 
attacked. Second, the results presented are expected to comprise a limiting case of those that ap- 
ply to more complicated problems. And third, preliminary results suggest that many of the results 
presented here apply generally to systems that are far from integrable. This topic will be discussed 
in more detail elsewhere. 

It is also worth mentioning that the results presented here are expected to apply to other fields, 
as well. Twinkling starlight is a familiar example, but far more exotic systems with the same 
dynamical foundation exist. Two such examples are the gravitational bending of light passing 
through galaxy clusters |}3], and the quantum mechanical waves associated with low temperature, 
conduction electron transport through semiconductor materials |33|. Recent experiments of this 
latter system show the electron transport breaking up into coherent channels that follow bundles 
of classical rays that are only weakly unstable in spite of the semiconductor acting as a random 
medium |Mj . This is a manifestation of the remaining stable and nearly-stable rays mentioned in 
Sect. II in connection with the KAM theorem; such rays will show up again later in this section. 

Ray stability analysis relies on the stability matrix Q defined in Eq. (|l3| ) from which a great 
deal of information can be deduced. It describes the behavior of all the rays that remain within 
an infinitesimal neighborhood of a reference ray over the course of its propagation. As seen in 
Eq. (|l~2|), geometric amplitudes for a point source are controlled by the matrix element (721; for a 
plane wave initial condition ^22 is the relevant quantity. For pedagogical purposes, however, it is 
simpler to consider the trace of the matrix; for chaotic rays, its exponential rate of increase and 
gross statistical properties are the same as those of any of the Q matrix elements. Since Q is 
diagonalizable by a linear, similarity transformation 



A = LQL^ 1 =► * " , (35) 




its trace is an invariant equal to the sum of the eigenvalues. The last form applies to systems with a 
single degree of freedom because the determinant is unity. Three distinct cases may arise. The first 
is |Tr(Q)| < 2 which is linked to stable motion, and it is customary to denote A = exp(i9r). The 
second case is |Tr(Q)| = 2, and it is often called marginally stable because it is the boundary case 
between stable and unstable motion. The third case represents unstable motion, and is characterized 
by |Tr(Q)| > 2. There it is customary to denote A = ±exp(^r) where v is positive and real. 
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In systems that are far from integrable, the vast majority of rays falls into this last category. 
Their largest Lyapunov exponent takes on an alternative form to Eq. ((2^), which can be expressed 
as 

ln(|Tr(Q)|) 



vl = lim 



(36) 



In the limit, r — > oo, effectively all rays converge to a unique value for v^. At finite ranges, each ray 
has different properties, and it is quite helpful to introduce the concept of a finite range Lyapunov 
exponent (sometimes known as stability exponent). For unstable orbits, Tr(Q) = A + 1/A~A with 
little inaccuray giving 

(37) 

r 

if necessary the more precise definition v = cosh _1 [Tr((5)/2]/r could be used, but it makes little 
difference for the discussion. Thus, in the limit, v approaches vl for every ray. The relationship 
between these two exponents is rendered evident by introducing an ensemble of potentials U each 
of which generates a chaotic dynamical system of fixed i>l. The ensemble average of v converges 
to vl without the necessity of the r —* oo limit, except at very short ranges. 

A very important property of these relations is that even small or modest fluctuations in v 
produce immense fluctuations in Tr(Q) which imply similar fluctuations in qi\- Although, long 
propagation ranges in chaotic systems imply vast numbers of eigenrays, the wave field may be 
dominated by those far fewer terms that appear in the summation with very small q<i\ ■ 

An analytic expression for the root mean square exponential rate of increase of Tr(Q) has been 
derived using techniques relying on Markovian assumptions |il 



^(Q)rms ~ exp(z/r) . 



where 



V 



d 2 U(z;r-£) 



dz 2 



d 2 U(z;r) 



z = z 



dz 2 



1/3 



(38) 



(39) 



Z = Z Q 



and the brackets, < ... >, denote ensemble averaging over different realizations of the potential U. 
Numerical evaluation of Eq. (39) for model systems with some realism is usually necessary, but 
analytic results are available for simplified models such as Gaussian random single scale potentials. 



Interestingly, v' is greater than the Lyapunov exponent, independent of the range involved. 
In fact, with increasing range it rapidly approaches a constant. The important distinction lies in 
whether the ensemble averaging occurs before or after taking the natural logarithm. The fluctua- 
tions are strong enough that v' > vl (= (z/)) or alternatively 



ln(|Tr(Q)| 2 ) / ln|Tr(Q)| \ 
2r \ r I 



(40) 



for all appreciable r. 

It turns out that, for a homogeneous background and single scale Gaussian random medium, 
to an excellent approximation except in the far tails, the probability density of v is a Gaussian of 
mean and variance 



(41) 



This distribution has also been observed in numerical simulations based on more realistic ocean 
models. This topic will be discussed in more detail elsewhere. Numerical calculations also suggest 
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that the Gaussian statistics are obtained for each realization of the random medium just by choosing 
different initial conditions; i.e. an ensemble of media is not necessary. A Gaussian density for v 
implies a lognormal density for |Tr(Q)| whose parameters are fixed by the Gaussian density's 
mean and variance. A property of lognormal densities is that any power of the variable also 
has a lognormal density. Thus, |Tr(<5)| 7 has a lognormal density as well; 7 = —1/2 relates to the 
semiclassical approximation. However note that depending on the value of 7, the lognormal density 
may fail as an approximation less far out in the tails. 

In the limit of long range, the fluctuations in |Tr(Q)| grow without bound despite v approaching 
vl for every ray. Just as there are highly unstable rays, there are also rays which are stable or nearly 
stable. The lognormal density gives a prediction for what approximate proportion is left for a given 
propagation range. Asymptotically, the proportion of nearly stable rays, whose stability exponent 
is less than some small value u c , decreases exponentially with range as (ao/47rr) 1//2 exp(— r/ao) 
where a = 2(V - v L )j (v L - v c ) 2 . 

The connection to the statistical behavior of wave field intensities arising from the lognormal 
density of classical ray intensities is not yet understood. There are a number of subtleties. To 
begin with, for long ranges of propagation, the most naive picture of semiclassical theory leads to 
the expectation that the wavefield is made up of an extremely large number of extremely small 
contributions. On average, if diffusive growth in phase space is neglected, the growth of the number 
of eigenrays and the decay of a typical ray intensity (squared amplitude) should occur at the same 
exponential rate. If this is true, energy conservation considerations dictate that at long range 
the constituent ray arrivals have effectively random phases. In other words, a set of N random, 
uncorrelated numbers of scale N~ 1 / 2 maintains an order unity summation as N — > 00. Chaotic 
dynamics is deterministic, and at best, the semiclassical phases generated by the dynamics are 
pseudo-random at long range. For shorter ranges, correlations amongst the magnitudes and phases 
could alter the statistical predictions. These correlations remain to be studied. Furthermore, the 
possibility mentioned earlier that one or few very strong terms at short range could dominate the 
sum increases the difficulty in finding a unique approach to an asymptotic statistical limit. 

A second difficulty is that the lognormal distribution has long tails indicative of its broad range 
of fluctuations. It has the form y~ a o^ n v which can be verified to approach zero for large or small y 
faster than any power of y, but does not decay exponentially. A sum of random numbers chosen 
from long tailed densities may not approach a standard central theorem limit (the Lorentzian 
density remains Lorentzian under repeated convolution), or may approach it very slowly. For the 
ray chaos problem, as more and more eigenrays exist with increasing range, the breadth of the 
density is also increasing. If the approach to a central limit theorem is slow enough, then it might 
never be reached under these circumstances; the limiting density still needs to be worked out. 

A third, extremely important difficulty is the appearance of caustics. Their number proliferates 
at the same exponential rate for chaotic systems as the number of eigenrays. At the caustics, the 
semiclassical expressions diverge and introduce infinities. They correspond to the extremely small 
values of instabilities in the tail of the lognormal density where breakdown of the statistical laws 
are likely. Thus, the lognormal expression does not, for example, capture the physics of singularity 
dominated fluctuations that are characterized by twinkling exponents which depend on the types 



of caustics present [36, g7|, 38]. Some incorporation of deviations from lognormality appears to 



be inescapable. In addition, the presence of exponentially proliferating numbers of caustics calls 
into question the very relevance of semiclassical methods and their usefulness in predicting the 
statistical properties of the wave field. This consideration is intimately linked with discussions of 
the validity of the semiclassical approximation for chaotic systems of which we give an overview in 
the next section. 

A final noteworthy complication arises in the analysis of sound fields produced by a broadband 
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source. Interference must still be accounted for, but only among ray arrivals whose travel times 
fall within intervals whose duration is the recriprocal bandwidth (A/) -1 of the source. This phe- 
nomenon is further complicated when pulse shapes and phase shifts at caustics - corresponding 
to Hilbert transforms in the time domain - are taken into account. Surprisingly, in spite of the 
widespread use of broadband sources in ocean acoustic experiments, we are unaware of any attempt 
to systematically account for all of these complexities in a theory of broadband intensity statistics. 



IV. WAVE CHAOS 

'Wave chaos' is the study of wave systems that, in the ray limit, exhibit unstable dynamics (i.e. ex- 
ponential divergence of neighboring rays); therefore the ocean acoustics problem can be thought of 
as a wave chaos problem. A completely analogous definition of the phrase 'quantum chaos' is also 
widely used. These two subjects are similar enough that much, if not most, of the progress in either 
domain carries readily over into the other; we therefore do not bother to distinguish them here. It 
turns out that there are significant conceptual difficulties with attempting to associate wave chaos 
with an unbounded, exponential growth of wavefield complexity in the hope that a straightforward 
generalization of the underlying ray chaos manifests itself. The shorthand explanation of this dif- 
ficulty is often crudely stated something like, "the finite wavevector (nonzero Planck's constant H) 
smooths over the intricate, fine scale details of the chaotic dynamics, not allowing them to be seen." 
The correspondence principle for chaotic systems, in fact, is quite subtle, and we avoid going into 
this fascinating subject except for the issue of the breakdown of semiclassical theories which we 
summarize next [39]. For further details though, we refer the reader to Ref. |4(J for a discussion 



of some aspects of wave chaos in underwater acoustics, and to some recent literature relating to 



quantum chaos [41, 42 1. 



The aspect of the breakdown of semiclassical theories which interests us most for the purposes 
of this paper is whether or not it is possible to construct, on solid theoretical ground, a ray-based 
theory valid beyond where ray chaos has fully developed. If not, ray-based predictions under chaotic 
conditions that match data in long-range propagation experiments would have to be considered 
accidental, and not an explanation of the essential physics of the problem. It may turn out that 
the eventual answer to this question is not unique, and depends upon which quantity one wishes 
to explain. For example, statistical predictions may be more robust than detailed, point-wise, 
wavefield predictions. On the other hand, if the breakdown occurs on a scale much longer than 
that of the development of ray chaos, there remains a great deal of theoretical work to be pursued. 

There is some hope for optimism in this regard. We begin by distinguishing between two classes 



of dynamical systems. The first class is that of simple, chaotic systems. Equations (28) define a 
system of this type. The equations of motion may be deceptively simple to write down, yet the 
solutions highly complicated, and, for all intents and purposes, analytically intractable (chaotic). 
The second class has complicated equations in the sense that the medium satisfies many of the 
criteria of randomness even though it is taken here to be deterministic (we may not know which 
determinisitic realization is given in a particular case) . One cannot take for granted the equivalence 
of the properties of these two classes of systems (chaotic versus random media), but we note that 
in certain limits there exist a number of common results (such as exponential ray instability with 
respect to initial conditions). 

In simple chaotic systems, about which more is known concerning semiclassical breakdown, we 
have to be careful to distinguish three levels of dynamics: classical, quantum, and semiclassical. We 
emphasize that the latter should be distinguished from the classical in that it takes classical infor- 
mation as input, but it actually generates an approximate construction of the quantum dynamics 
at the level of wave amplitudes and phases. The distinctions between time-evolving classical and 
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quantum expectation values of operators or classical and quantum probability densities have been 
studied since the development of quantum mechanics. Those quantities that correspond to each 
other initially are known to propagate similarly before diverging over a finite time scale called the 
Ehrenfest time. For simplicity, we shall focus only on the fact that one cannot delay the onset of 
interference phenomena in the quantum dynamics beyond the Ehrenfest time, and interference is 
necessarily excluded from the classical dynamics. For chaotic systems, the Ehrenfest time depends 



logarithmically on h [43], and beyond this time scale one finds all manner of complications such as 
an exponential proliferation of rays, increasing uncertainty whether the rays can even be calculated, 
and proliferating caustics in semiclassical theories. There is no debate whether the classical and 
quantum dynamics diverge beyond the Ehrenfest time (they do by definition); the relevant debate 
centers upon whether a ray-based, semiclassical theory can surmount these difficulties. 

There are various 'flavors' of semiclassical approximations possible. For example, the Wigner- 
Weyl calculus or other constructions of pseudo-phase spaces are ideally suited for exhibiting 
how the classical limit emerges from quantum mechanics in a semiclassical limit, but they are 
poorly adapted for describing interference phenomena. Indeed, mathematical proofs exist that 
such phase-space-based semiclassical approximations cannot be extended beyond a logarithmic time 
scale; see [^]. In this scenario, it can be fruitful to examine smoothed features of the wave field (as 
opposed to pointwise comparisons), such as can be obtained through the ray-based construction 
of the coarse-grained Wigner function [f5|. However, the development of semiclassical approaches 
that can be roughly described as time-dependent WKB theory handle interference naturally though 
multiple stationary phase (saddle point) contributions. This approach has been applied to a number 
of paradigms of chaos (the bakers map, the stadium, and the kicked rotor) with excellent, numerical 



results extending well beyond the logarithmic time scale [47, |48|, [49|] . Arguments leading to the 
expectation that the breakdown time depends algebraically on fi have also been presented [48, 49, 

Semiclassical breakdown in an idealized, but highly structured, ocean model consisting of a 
single realization of a ensemble with prescribed statistics has been investigated by Brown and 
Wolfson. They constructed the full semiclassical approximation using the classical dynamics and 
a Maslov-Chapman uniformization procedure. In their comparisons, they found the semiclassical 



approximation appeared to be working quite well beyond the onset of ray chaos [51]. This work and 
those quoted for simple chaotic systems are claiming that it is possible to extend a semiclassical 
theory for both chaotic and random media problems beyond a logarithmic time (range) scale, and 
thus, that ray-based semiclassical theories are viable candidates purporting to explain the essential 
physics of wave chaotic dynamical problems. 

V. CHAOS AND MODE COUPLING 

The connection between ray and modal expansions of acoustic wavefields in range-independent 



environments is well understood (see, e.g., Refs. [52|- |54{j ) and some generalizations have been 
derived |)^]-|57j for range-dependent environments. Thus, it should come as no surprise that there 
is a quantifiable connection between ray chaos and the modal description of the wavefield. A 



detailed description of this connection is described in Refs. pM] and [59|. In this section, the 



connection between ray chaos and mode coupling is briefly described. To simplify the presentation 
somewhat, we make use of a WKB analysis of the parabolic wave equation. Consistent with our 
use of the parabolic wave equation, the variables p' = cop, H' = cqHpe, U' = cqU and I' = cqI are 
used in this section after dropping the primes. Also, S = c^T is used in place of T. 

The principal result of this section is an approximate analytical expression for mode amplitudes 
in term of ray-like quantities. We discuss how typical features of ray chaos, such as exponential 
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growth (with range) of eigenrays contributing to the field point and coexistence of chaotic and 
regular ray trajectories manifest themselves in the mode amplitude range-dependence. The phe- 
nomenon of nonlinear ray-medium resonance that plays a crucial role in the emergence of ray chaos 
is shown to have an analog for modes which we call the mode-medium resonance According 



to the heuristic criterion proposed by Chirikov [19]-[21], chaos is a result of an overlap of different 
resonances. In terms of normal modes, chaos is shown to result from overlapping mode-medium 
resonances leading to complicated and irregular range variations of the modal structure. 

The normal mode representation of the solution to the parabolic wave equation ( |T7| ) is obtained 
by expanding ^f(z,r) in a sum of eigenfunctions of the unperturbed Sturm-Liouville eigenvalue 
problem |K| |6l| , 

1 d 2 ip m 

2 dz 2 

where ip m and E m are the eigenfunctions and eigenvalues, respectively: 



r. 



(42) 



(43) 



Here it is assumed that U(z,r) = Uq{z) + eV(z,r) and V is a measure of the magnitude of V. 
Each term in the sum (f43|) represents a contribution from an individual normal mode. In order to 
get simple semiclassical expressions for the amplitudes, 5 m , we project the ray representation of 
*&(z, r) of the form of Eq. (||) onto the normal modes. Assuming the latter to be normalized in 
such a way that 

/•oo 

dzipm(z)ip n (z) = 8 mn , (44) 



we need to evaluate the integrals 



Br, 



dz^(z,r)i) m (z) 



(45) 



Since we consider the geometric approximation to ^(z, r), it is natural to use the same approxima- 
tion for i/j m (z)- The corresponding formulas are usually referred to as the WKB approximations to 
the eigenfunctions [ 54 , 6C , . 

We shall assume that the potential Uq(z) is smooth, has only one minimum and its walls tend 
to infinity as z — > ±oo. When this assumption is combined with use of the WKB approximation, 
the eigenvalues of the action variable I m are determined by the quantization rule 



k Iri 



1 

m+~. 



(46) 



Then the eigenvalues of the "energy" are given by the relation E m = E(I m ), where the function 



E(I) is determined by Eq. (21). In the same approximation the m-th eigenfunction ip m ( z ) between 



its turning points can be represented as [54, 61] 



with 



i>m ( z ) = Qm(z) exp 



Qm(z) 



±i[ko dz'j2[E m -U (z')} 



7T 



fR^[2(E m -U (z)} 



1/4 ' 



(47) 

(48) 
(49) 
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where R m is the cycle length of the ray in the unperturbed waveguide with E = E m . In this form 
the eigenfunction is represented as a sum of two terms with rapidly varying phases and slowly 



varying amplitudes. Substituting Eqs. (47) and the ray representation of ty(z,r) of the form (g) 
into Eq. (45) yields a sum of integrals which can be approximately evaluated by applying the 
stationary phase technique [H] . This has been done in Ref. |p8| . Earlier, a related result was 
obtained in Ref. (6^]. Here we present only the final expressions for B m . 

It turns out that the amplitude of the m-th mode at the given range r is formed by contributions 
from the rays with action variables equal to that of the given mode, i.e., with 

I = Im- (50) 

Equation ( |50| ) singles out a set of rays which we shall call the eigenrays for the m-th mode. The 
action variable on the left is considered as a function of range, and initial values of the momentum 
Po and the coordinate Zq, i.e., 

I = I(p ,z ,r). (51) 

This function is determined by the solutions to Eqs. (||) and (|l9|), and Eq. (|2l|). For a point source 
the value of zq is the same for all rays; substitution of Eq.(|5l|) into Eq.(|50|) then yields an equation 
for pq, whose solutions define the initial momenta of the eigenrays contributing to the m-th mode. 

The process of identifying modal eigenrays is illustrated in the right panel of Fig. 2. Two curves 
are plotted: a surface of constant I and a segment of a Lagrangian manifold. The Lagrangian 
manifold corresponds, at r = 0, to a small angular aperture point source. Each intersection of the 
two curves corresponds to a modal eigenray. Eleven such intersections are shown. Under chaotic 
conditions this number can grow exponentially in range. (In contrast, in all range-independent 
environments there are two modal eigenrays for each mode, independent of range, whose launch 
angles have opposite sign.) 

The mode amplitude is given by the sum 

B m = J2 1 1 e (52) 

n yJko\dI/d Po \ po=PQn 

where each term represents a contribution from an eigenray with an initial momentum pQ n . The 
explicit expressions for the phase terms are (since the subsequent formulas describe characteristics 
of a single eigenray, the subscript n is omitted) 

<S> = S-S (z,I m )sgn(p) (53) 

and 

-*■*•> < M » 

Here z is the depth and p is the momentum of the eigenray at a range r, S and fi are its eikonal 
and Maslov index, respectively. 



Equations ([5l|)-(|54]) provide the analytical description of mode amplitudes in a range-dependent 
environment through the parameters of ray trajectories, i.e. through solutions to the Hamilton 
equations (0). These equations reduce the mode amplitude evaluation to a procedure quite analo- 
gous to that generally used when evaluating the field amplitude at the given point. This involves 
solving the Hamilton (ray) equations, finding the eigenrays, and calculating ray eikonals and some 
derivatives with respect to initial values of ray parameters. An important point should be stressed. 
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Although we expand the wavefield over eigenfimctions of the unperturbed waveguide, smallness of 
the perturbation has not been assumed. Our small parameter is the acoustic wavelength, 2ir/k, 
that should be substantially smaller than any physical scale in the problem. 

Having the comparatively simple expressions relating the mode amplitudes to rays, we can now 
discuss how the complicated ray trajectory dynamics reveals itself in the mode amplitude variations. 
For simplicity, we restrict our attention to a waveguide with a weak (e is now considered as a small 
parameter) periodic range dependence with a spatial period 2ir/Q, and we consider the case when 
only one mode is excited at r = 0, i.e. 

u(0, z)=il> m (z). (55) 



Analysis of the ray structure for this type of distributed source (see Refs. [58, ^]) shows that for all 
rays initial values of the action variable I are equal to I m . A situation which we call mode-medium 
resonance occurs when the value of I m satisfies Eq. (p4[). Due to the resonance, at ranges of order of 
1/Alv there will appear a bundle of rays with the action variables / in the interval \I — Iq\ < A/ max . 
This means (see Eqs.(|50|)) that starting with such ranges, the m-th mode is split into a group of 
2M modes with 



M = AI wax /k = 2y/eV/\uS\/k . (56) 

This expression is the modal analog of the standard expression (Eq. ^) for the width of a ray 
resonance. In the case of overlapping modal resonances it is natural to expect a further broadening 
of a group of modes. Moreover, as we discussed earlier, the overlapping of resonances causes the 
emergence of ray chaos with exponential proliferation of eigenrays. Under chaotic conditions the 
number of eigenrays contributing to a given mode also grows exponentially with range, giving rise 
to a very complicated range dependence of mode amplitudes. Numerical simulations presented in 



Refs. [58|, [5£]] support these statements. 

It might be assumed that that exponential proliferation of eigenrays contributing to a given 
mode leads to statistical independence of mode amplitude fluctuations under chaotic conditions. 
We expect, however, that the problem of describing mode amplitudes is considerably more rich and 
complicated. First, it should be recalled that generically the phase space of a chaotic Hamiltonian 
system contains both chaotic regions and "stable islands" formed by regular periodic trajectories. 
Some such regular rays will be eigenrays for some modes. Their contributions to modes cannot 
be considered as stochastic. Thus, we expect that under chaotic conditions there will be modes 
with amplitudes composed of two constituents: a chaotic one and a regular one. Numerical results 



illustrating this statement have been presented in Ref. 5£]. Another important phenomenon typical 
of chaotic dynamics, which may affect modal structure variations is stickiness, i.e., the presence of 
segments of a chaotic trajectory which exhibit almost regular behavior. The interval over which 
apparently regular behavior is observed can be long. In principle, one can presume that stickiness 
may cause some long-lasting correlations of mode amplitudes ffOfl . 

Our ray-based description of normal mode amplitudes has restrictions that are very much like 
those of standard ray theory. In particular, at some points the contribution from an eigenray 
to a given mode can be infinite. This occurs when the derivative in the denominator in Eq. 
( |52"|) vanishes. (The same comment applies to Eq. (|6C|), below.) Such divergences represent 
analogs of standard ray theory caustics. Under conditions of ray chaos the number of such caustics 
grows exponentially with range, spoiling applicability of the ray-based description already at short 
distances. This issue is discussed in Ref. |59[| . On the other hand, in Ref. |5J| (see also Ref. |46| 1) 
it was demonstrated numerically that the approach considered in this section can properly predict 
squared mode amplitudes smoothed over the mode number at ranges of order of at least ten inverse 
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Lyapunov exponents. This result is rather encouraging. It gives us hope that energy redistribution 
between modes can be comparatively easily analyzed using simple ray calculations at ranges of the 
order of a few thousand km. 

So far, we have considered only a cw field. For a signal with a finite bandwidth, the mode sum 
must include an integration over frequency a. The pulse signal at the point (z, r) can be represented 
as 

u(z, r,t) =Y^ u m (z, r, t), (57) 

m 

where 

u m (z,r,t) = das(a)B m (r,a)ip m (z,a)e V c o J (58) 



with s(er) being the spectrum of the initially radiated pulse. In Eq. ( |5"8| ) we indicate explicitly the 
dependencies of B m and ipm on a which have been omitted until now. Each term, u m (z,r,t), in 



the sum (57) can be interpreted as a pulse carried by an individual mode and we shall call it the 
"mode" pulse. The mode pulse, in turn, can be regarded as a superposition of elementary pulses 
representing contributions from different terms in the sum (j52|): 



u m (z,r,t) =\ das(a)G m (z,r,a)e \ c ° /, (59) 



i'n+r 



where 



2 cos (ko H dz' ^2[E m - U (z')\ - tt/4) 



G m {z,r,a) = = — — -j- . (60) 

/k r\dl/dp \ m=mn ^[2{E m -U (z)] 1 l i 



The above expressions depend on mode parameters with the subscript m and eigenrays parameters 
with the subscript n. Note that both types of parameters depend on frequency. Although the 
trajectories of eigenrays contributing to the given mode obey frequency-independent Hamilton 
equations (||), their starting momenta determined by Eq. ( |50| ) depend on the eigenvalue I m . But 
the latter, according to Eq. (|46| ), does depend on frequency. Under chaotic conditions the number 



of terms in Eq. (5£) can be huge leading to a very complicated shape of the mode pulse. 

The expectation that mode pulses are very complicated under chaotic conditions is consistent 
with the numerical results in Ref. IjSq ], where broadband parabolic-equation simulations of sound 
transmission through a deep water acoustic waveguide with inhomogeneities induced by random 
internal waves are described. It was shown that, due to mode coupling, mode pulses were several 
times longer than was the case in the background range-independent environment, and acquired 
irregular shapes. In contrast, from the ray perspective, the same weak inhomogeneities caused 
steep eigenrays to split into clusters of eigenrays (micromultipaths) whose travel time spreads were 
small and whose centroid had a travel time that was close to that of the eigenray in the background 
environment. (Recall also section II of the present paper.) The authors of Ref. |33| concluded that 
"while the high modes may be strongly affected by internal waves they are coherent enough that 
when they are synthesized together localized wave front results." A qualitative explanation of this 
phenomenon has been offered in Ref. [64]. In that paper it was shown that mode pulses may be 



considerably distorted due to mode coupling already at ranges so short that chaotic ray dynamics 
has not yet had a chance to reveal itself and every mode is formed by contributions from only two 
eigenrays. It was also demonstrated how distorted pulses carried by individual modes can combine 
to produce much less distorted ray-like pulses at the receiver. 

Although the mode coupling relations presented above are not easy to test experimentally, we 
emphasize that these results are important because of the insight they provide into the underlying 
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propagation physics. In this regard, it should be noted that the mode coupling relations presented 
above directly address the connection between ray chaos and finite frequency propagation effects, 
i.e., the wave chaos problem. 

VI. DISCUSSION 

In this paper we have reviewed results relating to ray dynamics in ocean acoustics. All of these 
results are intimately linked to the Hamiltonian structure of the ray equations. Most previous 
studies have emphasized the applicability of KAM theory to oceans with periodic range-dependence 
and the extreme sensitivity of individual chaotic ray trajectories. To complement these ideas, 
considerable attention has been focussed on oceans with nonperiodic range-dependence, and we 
have discussed some forms of stability of distributions of chaotic rays. Also, we have discussed 
subjects such as nondegeneracy violation, ray intensity statistics and the connection between ray 
chaos and mode coupling that either haven't, or have only recently, been explored in an underwater 
acoustic context. 

Although our understanding of ray dynamics is currently incomplete, it should be clear that 
the most pressing problem in this context is to better understand the connection between the 
ray dynamics and the corresponding finite frequency wavefields. For example, even the seemingly 
simple task of translating ray intensity statistics to wavefield intensity statistics is complicated by 
the necessity of making an additional assumption about relative phases and correcting for diffractive 
effects. 

It is our hope that the theoretical results presented in this paper provide a foundation for the 
analysis of measurements of sound fields at long range in the deep ocean. An analysis of this type 
will be presented in a forthcoming paper. Our long-term goal of developing tools that can be used 
to assist in the analysis of measurements accounts, in large part, for the considerable attention that 
we have devoted to ocean structures with nonperiodic range-dependence. 

The ideas and results that we have discussed differ in some important respects from more 
commonly applied ideas and results associated with the study of wave propagation in random media 
(WPRM). For example, most long-range underwater acoustic propagation WPRM theories (see, 
e.g., |6q1 ) assume that stochasticity is caused exclusively by internal waves. From the deterministic 
chaos point of view, this assumption is difficult to justify inasmuch as non-internal-wave (e.g., 
mesoscale) structure may excite ray-medium resonances and chaos. Also, in the deterministic 
chaos point of view, the excitation of ray-medium resonances generally leads to a mixed phase 
space. Most approaches to WPRM, on the other hand, invoke the assumption that the perturbation 
to the background sound speed structure is delta-correlated; this leads to stochasticity, but not to 
a mixed phase space. These topics will be discussed in more detail in future work. 

Finally, we wish to remark that the importance of ray methods is not diminished by recent 
advances, both theoretical and computational, in the development of full wave models, such as 
those based on parabolic equations. The latter are indispensible computational tools in many 
applications. In contrast, the principal virtue of ray methods is that they provide insight into the 
underlying wave physics that is difficult, if not impossible, to obtain by any other means. For this 
reason ray methods remain important. 
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